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1. Introduction. 

1.1. Polynomial transforms. Let b = {p (i),...,|)„_ 1 (x)} c C[x] be a lisfl 
of complex polynomials that form a basis of the space of polynomials of degree less 
than n, and let a = {ao, ■ ■ • , a„-i} C C be a list of distinct complex sample points. 
A polynomial transform is the matrix-vector product 7\ Q a;, where x € C™ and 7\ a 
is the n x n matrix whose (k,£)-th element is defined as pe(ak), < k,£ < n: 



V b , a = 



( Po(ao) 
Po(ai) 



Pi(ao) 
Pi(ai) 



Pn 
Pn 



-i(«o) \ 
-i(ai) 



\p (a n -i) pi(a n -i) 



(1.1) 



Pn-l(a n -l)J 



By a slight abuse of notation, we also refer to 7\ Q as a polynomial transform. 

Polynomial transforms are known in the literature under different names. For 
example, in [T3] and [5S], the authors refer to Vb, a as a discrete polynomial transform. 
In [20] . the authors call it a polynomial Vandermonde matrix. The most well-known 
example of a polynomial transform is the discrete Fourier transform (DFT). 

Polynomial transforms have a number of important applications. For example, 
they are used for interpolation and approximation |17j . solving differential equa- 
tions [6J, data compression and image processing pTrt24) . and the DFT specifically is 
widely used for spectral analysis and fast computation of correlation and convolution. 

The origin and main motivation for our work lies in the algebraic signal process- 
ing theory [29ll311f55] . This theory identifies polynomial transforms as equivalent to 



"This work was supported in part by NSF grant CCF-0634967. 

t Department of Electrical and Computer Engineering, Carnegie Mellon University, Pittsburgh, 
PA 15213 (asandryh@andrew.cmu.edu, pueschel@ece.cmu.edu). 

' Departments of Biomedical Engineering and Electrical and Computer Engineering, Carnegie 
Mellon University, Pittsburgh, PA 15213 (jelenak@cmu.edu). 

1 Hcrcaftcr, we view lists as ordered sets, i.e., without duplicate elements. 



2 



A. SANDRYHAILA, J. KOVACEVIC, AND M. PUSCHEL 



(generalized) Fourier transforms for shift-invariant 1-D signal models, and establishes 
a connection between these transforms and the representation theory of polynomial 
algebras. This connection has been used to algebraically derive many known and 
new fast algorithms for the DFT and discrete cosine and sine transforms (DCTs and 
DSTs) [3ll[32] extending early ideas by Nussbaumer [27] ■ All these algorithms are 
derived and represented as factorizations of the transform matrix into a product of 
structured matrices with low computational costs. 

In this paper, we develop a new algebraic method for a polynomial transform 
factorization. It is based on viewing the associated polynomial algebra as a regular 
module and decomposing it into an induction using a chosen subalgebra. This decom- 
position, performed in steps, yields a factorization of the polynomial transform. If all 
factors have sufficiently low computational costs, this factorization is a fast algorithm. 

Our method extends the approach in 31,32 to its most general form. As an 
application, we derive novel fast general-radix algorithms for the DFT and the DCT 
of type 4 that require only 0(n\ogn) operations instead of n 2 . 

1.2. Related Work. Over the last decades, decompositions that lead to fast 
algorithms have been studied for certain polynomial transforms. Among them, the 
DFT is arguably the most famous and well-studied. The discovery of the Cooley- 
Tukey fast Fourier transform (FFT) algorithm [TT] , which reduced the computational 
cost of DFT n to O(nlogn) operations, led to decades of research and numerous FFTs 
(see [3H1I10] and the references therein). 

Most other polynomial transforms of interest are related to the DFT and form the 
class of trigonometric transforms, since their entries are cosine and sine expressions. 
This class includes the DCT and the discrete sine transform (DST) of various types, 
as well as the real DFT and the discrete Hartley transform. Fast algorithms with 
O(nlogn) operations have been developed, for example in [3l[7lfT5 l l34 ll36ll42] . 

A more general class of polynomial transforms that were studied are those based 
on orthogonal polynomials [Hl[20l[28] . With the exception of DCT and DST, which 
belong to this group of transforms, the fast algorithms for this class reported in the 
literature require 0(n log 2 n) operation. 

Among hundreds of publications on this topic, most derived fast algorithms by 
clever, but often complicated manipulations of the matrix coefficients. This method 
provides little insight into the origin and the basic principles that account for the 
existence of these algorithms. 

Another thread of research that we refer to as an algebraic theory of transform 
algorithms has uncovered these principles for a large class of algorithms for trigono- 
metric transforms [30 , 32,41 . The theory exploits the connection between polynomial 
transforms and polynomial algebras and uses algebraic techniques to derive algo- 
rithms. As a result, most existing algorithms were identified as special cases of two 
basic theorems, the derivation is greatly simplified, and new algorithms were found. 

The origin of the algebraic approach is in [1,4,26,27], who recognized that the 
DFT„ can be interpreted as a decomposition matrix for the group algebra C[Z n ], 
where Z n is a cyclic group of order n QUE]- Since C[Z„] is identical to the polynomial 
algebra C[x]/(x n — 1), this decomposition is 

C[Z„] S C[x]/(x n - 1) C[x}/{x - uj n ) © • • • © C[x]/(x - C 1 ). (1.2) 

Algorithms are now derived by performing this decomposition in steps and reading 
off the respective matrices, which in turn factorize the DFT. 
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The group point of view was then generalized to derive fast Fourier transforms 
for group algebras C[G] for noncyclic finite groups G (SHHEHEIIISS] ■ Some of them 
were based on the induction for group algebras, a construction that is algebraically 
analogous to the method used in this paper. 

The polynomial algebra point of view was extended to derive and study larger 
classes of FFTs O^imEHEE] • 

The extension to the algorithm derivation of the full class of trigonometric trans- 
forms and a large class of algorithms was then accomplished in [321I4T] based on early 
ideas from [37][38]. Since all these algorithms are based on two theorems that gener- 
alize and account for the original Cooley-Tukey FFT, all the algorithms were called 
"Cooley-Tukey type." The close relation between transforms and algebra was fully 
developed and explained in the algebraic signal processing theory [3T]l3"3"]. 

In this paper we generalize the main theorem from |32j and hence the class of 
Cooley-Tuke type algorithms. Specifically, following the discussion in [32], we rigor- 
ously demonstrate in Chapter [5] that these algorithms can be viewed as based on a 
special case of algebraic induction. Then we generalize the construction method to 
its most general form and show that it produces novel algorithms. As examples, we 
derive new general-radix algorithms for the DFT and the DCT of type 4. 

2. Polynomial algebras and transforms. In this section we discuss polyno- 
mial algebras and demonstrate that their decomposition matrices are exactly poly- 
nomial transforms. We assume that the reader is familiar with the basic theory 
of algebras, modules, and matrix representations, even though we strive for a self- 
contained presentation in this paper. A good introduction to these topics can be 
found in [12,16,17 Below, we briefly review definitions and important properties. 

A vector space that is also a ring is called an algebra. In this paper, we work with 
polynomial algebras of the form A = C[x]/p(x). Elements of A are polynomials in x 
that are added and multiplied modulo p(x). We assume p(x) = YYk=o ( x ~ a k) S C[x] 
is a polynomial of degree n and separable, i.e. a k ^ ct m for k ^ m. A is a commutative 
algebra of dimension n with a multiplicative identity. 

A vector space M. that permits a multiplication by elements of A, such that 

am £ M. for any a € A, m e M., 

is called an A-module. The special case M = A is called a regular module. A 
subvector space Af < M. that is also closed under the multiplication by elements of 
A, is called an A-submodule of M. If M has only trivial submodules (i.e., {0} and 
itself), it is called irreducible. 

It follows from the Wedderburn theorem that a regular module M. = A can be 
decomposed into a direct sum of irreducible ^-modules [I2jll6]. This decomposition 
is accomplished by the Chinese Remainder Theorem: 

A: M -> ®V C[x]/(x-a k ), ^ 
s(x) i-> (s(a ) s(ai) ... s(a„_i)) . 

Suppose the basis of M. is a list of polynomials b — {pa(x), . . . ,p„_i(a;)}, and in 
each C[x]/(x — a k ) we choose the basis consisting of 1. Then the matrix that describes 
the isomorphism (|2.1|) is precisely the polynomial transform shown in (jl.lj) : 



Vb,a = Mafc)]o<M< 



(2.2) 
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Kind C Co, Ci C„(cos0) Symmetry Zeros (0 < k < n) 

1 st T l,x 
2 nd U l,2x 

3 rd V 1,2a; -1 

4 th W l,2x + l 

Table 2.1 

Chebyshev polynomials, their closed form C n (cos 9) , symmetry, and zeros. 



cos (n0) T_„ = T n cos (2fc 2 + 1)7r 

sinS — U "~ 2 COS n+1 

cos(n+^)9 „ _ (2fc + l)7T 

cos | - C ° h 2n + l 

sin(n+A)0 (2fe + 2)?r 

. a VK_„ = -VK„_i cos „T_ ' 



Namely, s(a;) = 53"=o s tPt{ x ) € becomes, in coordinate form, the column vector 

s(x) = (s si ... s„_i) T , 
and A(s(x)) in f|2 . 1 1) can be computed as the matrix-vector product 

A(*(a:)) = V b ,a ■ «(sc). (2.3) 



Example 2.1. If&={l,a;, . 



L } is the standard basis, then the polynomial 



transform (|2.2p is the Vandermonde matrix 

Vb, a = [ a k] < kt i <n ■ 

If, in addition, p(x) = x n — 1, then = w*, where w n 
the polynomial transform is precisely the discrete Fourier transform 

DFT n =[ W «] 



(2.4) 



with i = J — 1, and 



(2.5) 



Example 2.2. If 6 = {T (x|, . . . , T n _i(cc)} is the basis consisting of the Cheby- 
shev polynomials of the first kinqj, then the polynomial transform has the form 



Vb,a ~ [ T d a k)] < k j <r 

(2fc+l)7T 

2n 

mial transform is the discrete cosine transform of type 3 

(2fc+ 1)£tt" 



(2.6) 



If, in addition, p(x) = T n (x), then a k — cos ^^n™ ( see Table l2~Tj) . and the polyno- 



DCT-3„ = 



cos ■ 



2n 



(2.7) 



0<k,l<n 



Scaled polynomial transforms. The notion of a polynomial transform can 
be generalized by allowing a different choice of a basis in the C[x]/(x — a k ) in (|2.1[) . 
Namely, if we choose the basis {ck},c k <E C in each C[x]/(x — a^), then (|2.2[) becomes 
the scaled polynomial transform 



VI, 



diag — , 



Cn-l 



(2.8) 



2 Chebyshev polynomials Cfe are the polynomials that satisfy the two-term recurrence C k +i = 
2xC k — Ck—i |25| . Hence, the whole sequence of polynomials is determined by Co and C\. By 
setting x = cos 9, Chebyshev polynomials can also be expressed in their trigonometric closed form 
as functions of 9. These and other properties are shown in Table |2~T1 
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with Vb a as defined in 

Example 2.3. Letp(x) — T n (x), and choose the basis b = {Vo(x), . . . , V n -i(x)} 
in Ai, where Ve(x) is the £-th Chebyshev polynomial of the third kind. If we choose 
Cfc = 1/ cos ^ k+ 2n » then tne associated scaled polynomial transform is the discrete 
cosine transform of type 4: 



(fc+l/2)7T 

V b , a = diag < fc< „ cos 



cos • 



2n 

(k + l/2)(£+l/2)n 



(fc+l/2)(l+l/2) T 



cos^t 

2n 



0<kj<n 



0<kj<n 



DCT-4,, 



Note that all 16 types of discrete sine and cosine transforms are scaled or unsealed 
polynomial transforms with bases consisting of Chebyshev polynomials |31] . 

3. Subalgebra and its structure. In this section we discuss the structure of 
subalgebras of A = C[x]/p(x). 

3.1. Definition. Choose a polynomial r(x) € A, and consider the space of poly- 
nomials in r(x) with addition and multiplication performed modulo p(x): 

B={Y, c kr k {x) modp(x) | c k e c}, (3.1) 

k>0 

where all sums are finite. We call B the subalgebra of A generated by r(x) and write 
this as B = (r(x)) < A. 

3.2. Structure. Given r(x) € A, we first determine the dimension of B = 
(r(x)). Then we identify B with a polynomial algebra of the form C[y]/q(y) with 
a suitably chosen polynomial q(y). 

Let a = {ctQ, . . . , a„_i} be the list of roots of p(x). The generator r(x) maps a 
to the list P = {/?o, • ■ • , Pm-x\, such that for each oik € a there is a (ij € (3, for which 
r(ak) = (3j- Hence, m < n, since for some k and I we may have r[a.k) = r{ai). 

Theorem 3.1. The dimension of B = (r(x)) is dimS = m = \(3\. 

Proof. Let d = dim ,8. Since B < A, then dim£> < dim^l and the polynomials 
|l,r(x), . . . ,r ra-1 (cc)} span the entire B. From the isomorphism (I2.1[) we obtain 

d = rank (A(l), A(r(x)), . . . , A(r™ _1 (.T))) 
= rank [r £ K)] < M< „. 

Since r{otk) G P an d |/3| = ?n, the above matrix has only m different rows; hence, 
d < m. On the other hand, it contains the full-rank m x m Vandermonde matrix 

as a submatrix; hence, d> m. Thus, we conclude that d — dimS = m. □ 
Next, we identify B with a polynomial algebra. 

Theorem 3.2. The subalgebra B = (r(x)) can be identified with the polynomial 

algebra C[y]/q(y), where q(y) = Iljlo 1 (v ~ $i)> v ^ a ^ e following canonical isomor- 
phism of algebras: 

k: B -j. C[y]/q(y), 
r{x) i y y . 
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We indicate this canonical isomorphism as B = C[y]/q(y). 

Proof. Observe that B and C[y]/q(y) have the same dimension m, and n maps 
the generator r(x) of B to the generator y of C[y]/q(y). Hence, it suffices to show 
that q(r(x)) = mod p(x) in B. From (|2.1[) we obtain 

A(q(r(x))) = (q(r(a )) ... q{r(a n ^)) T 
= (0 ... o) T , 

which implies that q{r{x)) = mod in .4, and hence in B. □ 

Let c = {qo(y), ■ ■ ■ , qm~i(y)} be a basis of C[y]/q(y). The polynomial trans- 
form (|2.2j) that decomposes the regular module C[y]/q(y) (and hence the regular 
i3-module B) is given by (|2.1j) as 

^c,/3 = fe(^j)] <j/<m ■ 

Example 3.3. Consider the polynomial algebra A — C[x]/(x 4 — 1) with a — 
{l, — i, The polynomial r\{x) = x 2 generates the subalgebra B\ = (ri(x)) = 

C[y]/(y 2 — 1) of dimension 2, since ri(x) maps a to /3 = {l, — l}- 

The polynomial r2(x) = (x + x~ 1 )/2 — (x + x 3 )/2 generates the subalgebra 
£>2 = (r2(x)) = C[y]/(y 3 — y) of dimension 3, since r 2 (x) maps a to f3 = {l,0, — l}. 

4. Module induction. In this section we introduce the concept of module in- 
duction, which constructs an ,4-module M from a S-module N for a subalgebra 
B < A. We will show that every regular .A-module is an induction, which is the basis 
of our technique for polynomial transform decomposition. 

4.1. Induction. Similar to the coset decomposition in group theory [T21H6] . 
we can decompose a polynomial algebra A = C[x]/p(x) using a subalgebra B and 
associated transversal: 

Definition 4.1 (Transversal). Let B < A be a subalgebra of A. A transversal 
of B in A is a list of polynomials T = {to(x), . . . ,tL-i(x)} C A, such that, as vector 
spaces, 

L-1 

-4 = k[x)B = t {x)B © • • • 8 t L ^(x)B. (4.1) 

i=a 

Later, in Theorem l4.61 we establish necessary and sufficient conditions for a list of 
polynomials to be a transversal of B in A. In particular, for any B < A there always 
exists a transversal. 

Given a transversal of B in A, we define the module induction, which is analogous 
to the induction for group algebras in [12] . 

Definition 4.2 (Induction). Let B < A be a subalgebra of A with a transversalT 
as in (|4.1[) . and letM be a B -module. Then the following construction is an A-module: 

L-1 

M = 0^(x)7V, (4.2) 

£=0 

where the direct sum is again of vector spaces. Lt is called the induction of the B- 
module M with the transversal T to an A-module. We write this as M. — M A. 
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In this paper, we are primarily interested in regular modules. These are always 
inductions, as follows directly from (|4.ip and (|4.2I) : 

Lemma 4.3. Let B < A with a transversal T. Then the regular module A is an 
induction of the regular module B: 

A = BtrA. (4.3) 



4.2. Structure of cosets. We have established in (|3.2j) that the subalgebra 
B < A, generated by r(x) £ A, can be identified with a polynomial algebra C[y]/q(y). 
Next, we investigate the structure of each B- module ti(x)B in the induction (|4.3p . 

Consider a polynomial t(x) G A. As in Theorem 13.21 let r(x) map a to ft and 
let q(y) — rij=o (y ~ ft)- Further, let a' = {otk \ t(ctk) ^ 0} C a be the sublist of a 
that consists of those otk that are not roots of t(x). Finally, let r(x) map a' to ft C f3, 
and denote \/3'\ — m! . 

Theorem 4.4. The dimension of t{x)B is &\vc\t{x)B = |ft| = ml . 

Proof. The proof is similar to that of Theorem 13.11 The list of polynomials 
{t(x), t{x)r{x), . . . , t(x)r"~ 1 (x)} generates t(x)B as a vector space. Using the iso- 
morphism A in (|2.1[) we obtain 

dim (t(x)B) = rank (A(t(x)),A(t(x)r(x)), A(t(a;)r"- 1 (a;))) 

= rank[i(a fe y(a fe )] 0£M<n (4-4) 

= rank (diag (*(«fc)) 0$fc<n ' [ r 'K)] < w< „)' 

Theorem 13 . 21 shows that [r f (a^)] 0<fc f<n has exactly m = \f3\ linearly independent 
rows of the form 

(i /?, ^ ... /r 1 )- 

For each ft , the above row contributes exactly 1 to the rank of the matrix (|4.4[) 
if and only if there exists ctk such that f(afc) 7^ and r(oik) — ft- Since there are 
exactly |/3'| = ml such values of ft , we conclude that dim (t(x)B) — ml . □ 

Next, we identify the B-module t{x)B with a C[y]/g(y)-module. 

Theorem 4.5. The B-module t(x)B can be identified with the C[y] / q(y) -module 
C[y]/q'(y), where q'{y) — Y\.b eB' (v ~ ft)' v ^ a ^ e m °dule isomorphism 

v: t(x)B -+ C[y]/q'(y), 

t(x)r k (x) i-> y k . ^ ' ' 

By a slight abuse of notation, we write t(x)B = C[y]/q'(y). This is an isomorphism of 
modules and should not be confused with the isomorphism of algebras in Theorem \3.£[ 

Proof. It follows from Theorem 14.41 that H(x),t(x)r(x), . . . ,t(x)r m is a 

basis of t(x)B, viewed as a vector space. On the other hand, (1, y, . . . , y m ~ x } is 
obviously a basis of C[y]/q'(y), also viewed as a vector space. Hence, r] in (|4.5I) is a 
bijective linear mapping between t(x)B and C[y]/q'(y). 

In order for 77 to be an isomorphism of modules, it must also be a module 
homomorphism — it must preserve the addition and multiplication in t(x)B and C[y]/q'(y). 
Namely, for h(x) £ B and u(x),v(x) £ t(x)B, the following conditions must hold: 

n(u(x) + v(x)) = r)(u(x)) + rj(v(x)), 
T](h(x)v(xy\ = n(h(x)) ■ n(v(x)). 
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The first condition is trivial. To show that the second condition holds, let h(x) = 
Tl=o h k r k (x) e B and v(x) = 1 Vjt(x)r^(x) G t(x)B. Then 

m+m' — 2 j m-\-m' — 2 j 

r](h(x)v(x)) = r)( ^2 ^2h k Vj- k t(x)r J ' (x)) = ^ 
j=a fe=o j=a fc=o 

m— 1 m' — 1 

fc=0 3=0 

Hence, ij is a module isomorphism. □ 

Note that, depending on t(x), the dimension of t(x)B may be smaller than the 
dimension of B: m! < to. This effect is called annihilation. 

Also, the definition of r\ in (|4.5|) assumes the standard basis {l, y, . . . , y m _1 } 
in C[y]/V(t/). If another basis {bo(y), . . . , b m '-i(y)} were desired, the corresponding 
basis in t(x)B would be {t(x)bo(r(x)), . . . , t(x)b m i -\(r(x))\ . 

As a consequence of Theorem 14.51 and the above discussion, decomposing the 
£>-module t(x)B with basis {t(x)qo(r(x)), . . . , t(x)q m /-i(r(x))\ is the same as decom- 
posing the C[y]/g(y)-module <C[y\/q'(y) with basis c = {q (y), ■ ■ ■ , q m '-i(y)}- The 
decomposition matrix is the same as for the regular module C[y]/q'(y) with the same 
basis, namely 

'Pc,p> = [<lt(!3 j )]o<M< m >- ( 4 - 6 ) 

4.3. Existence of a transversal. Consider T = {to(x), . . . , tL-i(x)} C A, and 
let dim(U(x)B) = m e for < £ < L. Then {t t (x), t e (x)r{x), . . . , t i (x)r mt - 1 (x)} is a 
basis of t((x)B, as follows from Theorem 14.41 Hence, T satisfies (|4.1|) if and only if 
too + • • • + m,L-i = n and the concatenation of bases 

L-l 

b'=\J {t t (x),...M*)r mt -Hxj} (4-7) 

is a basis in A. The following theorem states this condition in a matrix form. 

Theorem 4.6. Using previous notation, T is a transversal if and only if the 
following is a full-rank n x n matrix: 

M' = (D B | D 1 B 1 | ... | Di_iB£_i) , (4.8) 
where D e = diag {tt(a k )) < k<n , and B t = [r* (a*)] < fc< „ i0 < j<m< ■ 



Proof. The proof is similar to the proofs of Theorems 13.11 and 14.41 Observe that 
the fc-th element of b 1 in (|4.7p is mapped to the /c-th column of M' in (|4.8|) by the 
isomorphism A in (|2.1j) . Hence, 6' is a basis in A if and only if M' has exactly n 
columns and rankAf = n. □ 

It follows from Theorem 14.61 that for any algebra A and its subalgcbra B there 
always exists a transversal. For example, we can choose T = {to(x), . . . 
where te(a k ) = for £ ^ k and ti{ai) ^ 0. In this case M' = diag (tt (ai )) 0< £ <n 
in (|4.8|) is a full-rank diagonal matrix. 

Example 4.7. Consider the subalgebras constructed in Example 13.31 

For B\ = (x 2 ) of dimension 2, we can choose the transversal T = {l,x}, since 
|l,a; 2 } U {x,a; 3 } is a basis for A. Since x maps a to {l, — i, — we have a' = 
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{l,-i,-l,i} and^' = {1,-1}. Hence, q'{y) = (y-l)(y+l) and xBi = C[y}/{y 2 ~ 
is of dimension 2. 

For B2 = ((x + x~ 1 )/2) of dimension 3, we can choose the transversal T 
{l, (x — x _1 )/2}, since the corresponding matrix 



1) 



M' 



( 1 
1 

1 -1 

V 1 



1 1 
1 



\ 



/ 



from (|4.8[) has full rank. Since (x — x 1 )/2 maps a to {0, — z,0,«}, we obtain a' = 
{ — f3' = {0}, and thus q'(y) — y. Hence, (x — x~ 1 )/2 ■ B2 — C[y]/y is of 
dimension 1. 

5. Decomposition of polynomial transforms using induction. In this sec- 
tion we use the induction (|4.3I) to express the polynomial transform of A via the 
polynomial transforms of each te(x)B = C[y]/q' g (y) in (|4.1I) . 

As before, we consider A = C[x]/p(x), where p(x) — Yik=o ( x ~ afc ) - We y i ew ^ 
as a regular ,4-module with the chosen basis b = {po(x), . . . ,p„_i(x)}. 

Let B = (r(x)) < A be a subalgebra generated by r(x) € A, and B = C[y]/q(y) 
according to Theorem 13. 2[ where q(y) — (u ~ Pj) an< ^ P — {A)> • • • , Pm-i} ■ 

Suppose T = {t {x), . . . ,tx,_i(x)} is a transversal of B in A. Let each ti(x)B 
in (|4.ip be identified with a C[y]/q(y)-module C[y]/q( e \y) according to Theorem 14. 51 
where gW(y) = (tf - ft) and m/ = The basis = {&<%), . . . , bjg^fc/)} 

of C[y]/g^(y) corresponds to the basis {^(a;)6^(r(x)), . . . , ^(x)6^ f _ 1 (r(a;))} oiti(x)B. 
Hence, the corresponding polynomial transform (|4.6[) is o(<) . 

Theorem 5.1. Given the induction (|4.3I) . £/ie polynomial transform Vb, a can be 
decomposed as 



L-1 



V b>a = (D M I D X M X 



Dl.xMl-i) (07> 6Wii9 (o).B. (5.1 



Here, B is the base change matrix from the basis b to the concatenation of bases 

L-1 



(J {^(^(Ka:)),...,^^.!^))}. 



i?ac/i = diag (te(ctk)) 0<k<n is a diagonal matrix. Each is an n x mi ma 
trix whose (k,j)-th element is 1 if r(afc) is 
otherwise. © denotes the direct sum of matrices 



(V, 



L-1 

£=0 



fe(«),/3(«) 



&<°),/3C0) 



io i/ie j-th element of (3^ , and 



b(i),/3<i 



Proof. We prove the theorem for L — 2; that is, for .4 = to(x)B ti(x)B. The 
proof for arbitrary L is analogous. 
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Let B = C[y]/q(y) according to Theorem 13.21 where q(y) = YYJLo 1 (v ~ 
and (3 = {fo,...,p m -i}. For t e {0,1}, let t e (x)B = C[y}/q^(y) according 
to Theorem where q^(y) = U/n-epW (V ~ Pj) and ™t = Let = 

{&<%), . . . , b^M} be a basis of C[y}/q^(y). 

Let ^(a;)6W(r(x)) = {t e (x)b £) (r(x)), . . . ,^(x)6^ ) o _ 1 (r(a;))}. As we established 
in Theorem 14.61 6' = ioW>^ ( r ( x )) U (K 2 ')) is a basis of A. The original 

basis b can be expressed in the new basis b' as Pk{x) = J^T^q 1 Bkjto(x)bf'\r(x)) + 
Sfco" 1 Ck,iti(x)b^~\r(x)) . Hence, if B is the base change matrix from b to b', then 

Pb,a = T>v, a ■ B. (5.2) 

The ^-th column of B is (B ,g, £ mo _i^, C /, . . . , C TOl _i^) . 
Next, observe that 

(a;)6(°)(r(a:)),a I ' ti(x)b( 1 '> (r(x)),a ) • (5-3) 

For each £, the (fc,j')-th element of 'Pt«(x)6( i )(r(x)),a i s 

t/(a fc )6W(r(a fc )). Hence, 

7 , t f (a:)fc(*)(r(2;)),a = D l ■ M(. ■ V b {t) ^(t) , (5.4) 

where Af^ is an n x to^ matrix whose (k, j)-th element is 1 if r(a.k) equals to the j-th 
element of /3^, and otherwise; and Dp = diag [tAotk) ) 

V /0<fc<n-l 

Hence, from (|5.21I5.4[) we obtain the desired decomposition: 

Pb,a = {D M I D1M1) • (n(o)^(o) © 7> 6 <i),/jm) • B. (5.5) 

□ 

Corollary 5.2. Consider the n x m matrix M whose (k,j)-th element is 1 i/ 
r(afc) = fij and otherwise. Then 

1. M contains exactly n Is and n(m — 1) Os . 

2. Each matrix Mg in Theorem I5.il is a submatrix of M. It contains the j-th 
column of M if and only if (3j € ■ 

3. If the number of non-zero elements in the j-th column of M is Cj, then there 
are precisely Cj matrices among Mo, . . . ,Ml-i that contain this column. 

Discussion. The three factors in (15.11) correspond to the decomposition (12. lj) of 
the regular module A = A4 = C[x]/p(x) in three steps: 

Step 1. A is represented as an induction (|4.3p by changing the basis in A to the 
concatenation of bases b^ of tg(x)B, using the base change matrix B. 

Step 2. Each tg(x)B is decomposed into a direct sum of irreducible Z?-submodules, 
using the corresponding polynomial transform V^i) «(«). 

Step 3. The resulting direct sum of irreducible £>-modules is decomposed into a 
direct sum of irreducible .A-modules, using the matrix M. 

The factorization (|5.ip is a fast algorithm for 7\ Q if the matrices B and M have 
sufficiently low costs, since the recursive nature of the second step allows for repeated 
application of Theorem 15.11 We illustrate this with two examples of novel algorithms 
derived using this theorem in Section [6] 

Special case: factorization of p(x). A special case of Theorem 15.11 has been 
derived in [3171132) . Namely, assume that A — C[x]/p(x), and we can decompose p(x) = 
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q(r(x)). Then B = (r(x)) = C[y]/q(y), and any basis t = {l, ti(x), . . . , t k _i(x)} of 
C[x]/r(x) is a transversal of B in A. This leads to the following result. 

Corollary 5.3. Choose c = {co(y), . . . , c m _i(y)} as the basis ofC[y]/q(y). De- 
note the roots ofr(x)-/3 3 as 7W = {7^ , . . . , 7^} . Notice that U™^ 1 {7^ , . . . , 7^} 
is simply a permutation of {ao, . . . , a n -i} , and denote the corresponding permutation 
matrix as P. Then, the polynomial transform decomposition (|5.1|) has the form 



V b ., 



P- 



3=0 



L 



(5.6) 



Here, <E> denotes the tensor product of matrices. 

Corollary 15. 31 has been used to derive a large class of fast algorithms for real and 
complex DFTs, and DCTs and DSTs P PUSI] . Theorem O further generalizes this 
approach, and, as we show in the following example and in Section [SJ also yields fast 
algorithms not based on Corollarv l5.3l 

Example 5.4. Consider the polynomial algebra A — C [x] / (x 4 — 1 ) with basis b = 



{1, 



1 } . As we showed in Example 12.11 the corresponding polynomial transform 



is V b , a = DFT 4 . 

We continue from Example 14.71 First, consider B\ = (a; 2 ) and the induction 
A = Bi © xB\. Let us choose 6 (0) = {l,y} as the basis of C[y]/(y 2 - 1) = B x ; it 
corresponds to the basis {l,x 2 } of B\. We then choose i>W — {l,2/} as the basis 
of C[y]/(y 2 — 1) = xB\\ it corresponds to the basis {a;, a; 3 } of xB\. According to 

Theorem 15. 11 D a = diag ^1, 1, 1, 1J, D\ = diag ^1, — i, — 1, ij, 



M = Mi = 



V 1 / 



6«», /3(°: 



= V, 



= DFT 2 , 



and B is the base change matrix from {l, 2, a; 2 , a; 3 } to {l,x 2 } U {a;, a; 3 }. Hence, 



DFT 4 = 



-1 



V 



DFTg 



DFT, 



V 



(5.7) 



1 / 



As we show in Section MM (E3 is exactly the Cooley-Tukey FFT for DFT 4 [IT] . 

Next, consider £? 2 = ((x + a;^ 1 )/2) and the induction .4 = S 2 ffl(a;— a; _1 )/2-02' Let 
us choose &(°) = {T (v),ri(i/),T 2 (tf)} = {l,y,2y 2 - 1} as the basis of C[y]/(y 3 -y) S 



$2; it corresponds to the basis (1, (x + x 1 )/2,(x 2 



! )/2} of B 2 - We then choose 



M 1 ) = {l} as the basis of C[y]/y = (x — x 1 )/2 ■ B2] it corresponds to the basis 
{(a; - a;- 1 )/2} of (x - x- x )j2 ■ B 2 . According to TheoremEU D = diag (l, 1, 1, 1 

£>i = diag (0, -i, 0, i) , V bW pw = (1) = DST-li , 



/ 1 



M = 



\ 



,Mi 



1 

V 1 / 



M°),/3(0) 




DCT-I3, 
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and B is the base change matrix from {l, x, x 2 , a; 3 } to {l, (x + x 1 )/2, (x 2 + x 2 )/2}u 
{(x - x- 1 )/2}. Hence, 



DFT 4 



1 -i 
1 

V i i ) 



DCT-1 



3 



DST-li 



1 1 
1 

Vi -i ) 



(5.8) 



As we show in Section 16.41 (|5.8I) is the Britanak-Rao algorithm for DFT4 [5] . 

6. Fast Signal Transforms. In this section we apply the module induction 
to the construction of novel fast algorithms for trigonomatric transforms, which are 
the most important polynomial transforms used in signal processing. The efficient 
computation of these transforms is of crucial importance in most applications, and 
makes straightforward computation using 0(n 2 ) operations prohibitive. As mentioned 
in the introduction, many 0(n log n) algorithms have been derived for these transforms 
(e.g., [3 7 15,34 36,42 ) and the origin of these algorithms was revealed by the algebraic 
approach in |30U32[|4"T] . which also produced new algorithms. 

In this paper, we complete this work through Theorem 15.11 and its application. 
Specifically, we will derive two novel 0(n log n) general-radix algorithm that could 
not be obtained with the prior algebraic theory. 

We will first briefly touch on the algebraic signal processing theory to explain 
why these transforms are associated with polynomial algebras. Then we derive the 
Cooley-Tukey FFT as special case of Theorem 15.11 which motivates why we call all 
such algorithms "Cooley-Tukey type." Then we derive the novel algorithms, both 
of which generalize existing algorithms that had no satisfying algebraic explanation 
before. 



6.1. Algebraic Signal Processing. In 30 33J, the authors introduced an ax- 
iomatic approach to the signal processing called the algebraic signal processing theory. 
They observed that the basic assumptions used in signal processing are equivalent to 
viewing filters as elements of an algebra A, and signals as elements of an associated 
A- module M.. In particular, in the shift- invariant signal processing of finite discrete 
one-dimensional filters and signals A = Ai = C[x]/p(x) is necessarily a polynomial al- 
gebra. The choice of A = A4, together with a bijective mapping $ that maps samples 
from C™ to signals in A4, defines a signal model (A, A4,Q). 

The fundamental tool in signal processing is the Fourier transform, which com- 
putes the frequency content of a signal. From the algebraic point of view, the Fourier 
transform for a signal model (A,A4,&) is precisely the decomposition (|2.ip . It can 
be computed as a matrix-vector product (|2.3p with the appropriate polynomial trans- 
form (1231). 



6.2. Notation. Hereafter, we use the following special matrices: 
/„ is the identity matrix of size n. 

J n is the complimentary identity matrix of size n: its (fc, n — 1 — fc)-th element is 
1 for < k < ii, and otherwise. 

l n = (l 1 ... l) is a column vector of n ones. 
Z n is the n x n circular shift matrix: 

1 

I-n-l 
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where k divides n, is an n x n permutation matrix that selects elements of 
0, 1, . . . , n— 1 at the stride k; the corresponding permutation is ik+j i— > jm + i, where 
< i < m and < j < k. The (i, j)-th element of L k is 1 if j = | ^ fc (" +1 ) j m od n, 
and otherwise. 

K k = (I m © J m ffi I m © . . . )L£ , where k divides n, is another permutation matrix. 

TJ} — diag (w 1 ^ | < « < fc, < j < to^ , where the index i runs faster, and n = 

km, is a twiddle factor matrix used in the Cooley-Tukey FFT. 
Complimentary direct sum: 



v At, 



6.3. Cooley-Tukey FFT. We derive the general-radix Cooley-Tukey FFT us- 
ing Theorem 15.11 As was shown in [3D] , Corollary 15.31 is sufficient in this case. 

Consider A = M = C[x]/(x n - 1). Let b = {l,x, .. .,x n - x ) be the basis of M. 
As we showed in Example 12. 1[ the corresponding polynomial transform is DFT„. 
Assume n — km. Let r(x) — x k , and B = (r(x)). Then x l B = C[y]/(y m — 1), for 
I = . . . k - 1, and A = ®^x e B. Choosing the same basis b w = {1, y,... , y m_1 } in 
each C[y]/(y m - 1) = x l B yields V m ^m = DFT m . By TheoremEU we obtain 

DFTfcm = M • (4 ® DFT m ) • 5. 

Here, B = L k k m and M = (D M \ . . . \D h _ x M ) , where M = l k <g> I m , and Z) £ = 



' 0<j<km 

km ( j k> I \ L' r r \ rrikm j km 



diag ( wf- 7 ) for < ^ < fc. Hence, we can rewrite 

V Km /0<7<fcm 



M = Lf" (7 m ® DFTfc) T k L n 

to obtain the well-known general-radix Cooley-Tukey FFT algorithm [TTJ[32] : 

DFT fem = L k k m (I m ® DFT fc ) T fc fcm i^ m (J fe ® DFT m ) i£™ 

= L km (I m ® DFT fc ) T fc fcm (DFT m ® 7 fc ) . (6.1) 

6.4. New Fast Algorithms. In this section, we derive novel fast general-radix 
algorithms for DFT and DCT-4. Each of them requires 0{n log n) operations. To the 
best of our knowledge, these algorithms have not been reported in the literature. 

General-radix Britanak-Rao FFT. In [8], Britanak and Rao derived a fast 
algorithm for DFT2,„ that can be written as the factorization 

DFT 2m = X^{l m © Z- 1 )D* n (DCT-l rn+ i © DST-U-^B*". 

Matrices D^ 1 , B^ 1 , and X^ 1 are specified in (|A.4tiA.6[) by setting k = 1. 

In Appendix \K\ we derive the following general-radix version of this algorithm: 
Theorem 6.1. 



DFT 2fcm = Lf m (l 2m ® DFT fc )x^™L^™(/ m © Z^ 1 © /. 



2(fe-l)m 



2/r 



X 



DCT-Vi © DST-l m _i © I k -i ® (DCT-2 m © DST-2 m ) IB. 



2A-/, 
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Here, is a diagonal matrix, and B 2 ^™ 1 and X^ m are 2- sparse matrices (that is, 

with each row containing only two non-zero entries) specified in {A.$A.(fy . 

This factorization is obtained by inducing a subalgebra B — ((x k + x~ k )/2^ of 
an algebra A — C[x]/(x 2km — 1) with transversal to(x) = 1, ti(x) — (x k — x~ k )/2, 
t 2 j(x) = x j {x k + l)/2, and t 2j+ i(x) = x j {x k - l)/2 for 1 < j < k. 

DFT fc requires 0{k log k) operations; DCT-l m+1 , DST-l m _i, DCT-2 m , andDST-2 r . 
require 0(m log m) operations each |30[l32j. D 2 ^ 171 requires n — 2km operations and 
Bm m ana - -^m™ eacn require 3n operations. Hence, the algorithm for DFT„ in The- 
orem |nH] requires O(nlogn) operations. 

General-radix Wang algorithm for DCT-4. In [42], Wang derived a fast 
algorithm for DCT-4 2m that can be written as the factorization 



DCT-4 2m = K 2m ■ 



m_1 2m-2i'-l / ,\j r 2j + l-2m 

COS 5-7 7T ( — l)"'COS^ LI S 7T 

8m rt ^ / 8m 



cos 2j + 1 ~ 2m ^ (-l)J+ 1 coB 2m -^'- 1 7r 

8m v ' 8m 

2m\T _ | r 2(m-l) 



x(DCT-3 m ® h)(K z 2 m Y ■ Li (m - L> ■ / m _! (g) DFT 2 



1; 



In Appendix [B] we derive the following general-radix version of this algorithm: 
Theorem 6.2. 

DCT-4 2fem = K 2km (K 2m <g> DCT-4 fc )F„ 2fem • (DCT-3 m ® L 2fe )(X 2 " fe ) T 
X/*® | L^™" 1 ) • ® DFT 2 J (tf 2 2 £ m ) T . 



Here, Y 2km is a 2-sparse matrix specified in (IB.4[) . 

T/iis factorization is obtained by inducing a subalgebra B = (T 2 k(x)') of an algebra 
A = C[x]/T 2 km{x) with transversal t 2 j(x) = Vj(x) and t 2 j+i(x) = Wj(x)(V 2 k-i(x) — 
V 2k (x))/2 forO<j <k. 

DCT-4 fe requires O(fclogfc) operations, and DCT-3 m requires 0{m log m) opera- 
tions [501l32j . Y^ rn requires 3n operations, where n = 2km. Hence, the algorithm for 
DCT-4„ in Theorem 16 . 2 1 req uires O(nlogrt) operations. 



7. Conclusion. We have introduced a new approach to the factorization of poly- 
nomial transforms 7\ Q based on the decomposition of the underlying regular module 
A = M. = C[x]/p(x) into an induction. This approach is in its most general form since 
the underlying Theorem 15.11 allows for arbitrary subalgebras. Not every factorization 
based on this theorem yields a fast algorithm: it depends on the computational costs 
of matrices B and M that occur in its recursive application. 

However, we have shown that the theorem produces at least two novel general- 
radix algorithms for the DFT and a DCT. Both algorithms cannot be obtained using 
the prior Corollary 15.31 In addition, both generalize algorithms from the literature, 
which now become the special cases of radix 2. 

Equally important, we make another step towards a complete algebraic theory of 
fast algorithms for polynomial transforms. 

Future work. In addition to the DFT, DCT, and DST, other polynomial trans- 
forms have been studied. In particular, polynomial transforms based on orthogonal 
polynomials have found applications in such areas as function interpolation, data com- 
pression, and image processing 22-24 . For practical applications, fast algorithms for 
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this class of polynomial transforms are needed. With the exception of DCT and DST, 
the fastest algorithms, reported in the literature to date, require 0(n log 2 n) opera- 
tions (in particular, more than 43nlog2 n for n = 2 k ) [14,28 . The question is whether 
our approach can improve this bound for some or all of these transforms. 

Appendix A. Proof of Theorem 16.11 Consider A = M = C[x]/(x 2km — 



1), with basis {l, x, . . 



„2fcm-l 



} and afc = ^2km- The corresponding polynomial 



transform is DFT 2 fc m . 

By Theorem 13.21 the polynomial r(x) = (x k + x~ k )/2 generates the subalgebra 

B = (r(x))=C[y}/2(y 2 -l)U m ^(y). 

If we choose {Te(y)} 0<g<m+1 as the basis , the polynomial transform is 

DCT-Wi. 



2>(cos — ) 
to . 



0<k,£<m+l 

By Theorem 021 the S-module (x k - X - k )/2 ■ B = C[y]/U m -i{y). If we choose 
the basis {Ue(y)} 0<i<m _ 1 , then the polynomial transform is 



l^(cos — ) 
to 



0<k,Km- 



= diag ^1/s 



(k+ 1)tt 
sin — 



rn 



■ DST-1, 



DST-1' 



0<fc<m-2 



(A.l) 

Similarly, the ^-module x 3 (x k + l)/2-B C[y]/2{y-l)U m - X {y) for any 1 < j < k. 
If we choose the basis {V«(y)} c 



t/ / kit 
VWcos — ) 

TO 



0<k,£<m 



I 0<i<m' 

diag (l/ 



then the polynomial transform is 
kit \ 

"■OS 



2to/ 0<k<m 



• DCT-2 m = DCT-2' 



(A.2) 



Finally, the S-module x 3 (x k -l)/2-B = C[y]/2(y + l)[/ m _ 1 (j/) for any 1 < j < k. 



If we choose the basis {W^(y)} 



o<e<m 



then the polynomial transform is 



(fc+l)7T 

W^(cos- — ) 



0<k,£<m 



= diag ^1/ sin 



(fc+l)7T 
2 TO 



• DST-2 m = DST-2' 



0<fe<m 



(A.3) 

Using Theorem 14. 61 we can verify that to(x) = 1, t\(x) = (x k — x k )/2, t2j(x) = 
x J (x + l)/2, and t2j+i(x) = x 3 (x k — l)/2 for 1 < j < k, is a transversal of B in A. 
Hence, by Theorem 15. 11 we obtain the factorization 



DFT. 



M DCT-1 



m+l 



DST-l^.! 4-i ® (DCT-2^ © DST-2JJ )B, 



ylknl 



Here, B 2km is the base change matrix from {x e } Q ^,„^ , to the concatenation of bases 
of tj(x)B, < j < 2k, and by construction 



I 0<l<n-l 



Tj2km 
&™ — 



I,, 



J: 



\ 



m — 1 



l 



fe-i< 



/ 1 1 



-1 1 



V I'm — 1 Jm — 1 / 

M is constructed as follows. Let 

M = 



Jm-l 
Jrn — 1 



■ L 



2k i, 



(A.4) 



Im-l 



\ Jm-l / 
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Let Mq(Jo, . . . , jg) be the subset of columns of M with indices jo, . . . , jf, and let 

, for < j < 2k. Then 

0<i<n 

M = (D M I D X M X I D 2 M 2 I • • • I D 2k -iM 2k A , 



Dj = diag I tj (a 



where M x = M (l, . . . , m - 1); M 2 j = M (Q, . . . , m - 1) and M 2 j+i = M (l, . . . , m) 
for 1 < j < k. We can further rewrite M as 

M = Lf m {I 2m ® DFT fe )A^ m L^ m (/ m © Z" 1 © 7 a(fc _ 1)jn ). 

Here, matrix X^ m has the structure 



/4 



A m-l/i 



3 in — 1 



D: 







2m-l 



F 



(A.5) 



where 



q = l©diag(^ fem K n + l)/2) i 
D 3 = ((4m ~ w^)/2) © diag (^L>2 m - l)/2) i 



i<f<fc 



F = 1 



© diag ( - u4 fc ) i 



<i<fc 



After the substitution of DST-l' ra _ 1 , DCT-2^, and DST-2^ with DST-l m _i, 
DCT-2 m , and DST-2 m using (|A.HIA.3p . and simplification, we obtain the factorization 

DFT 2 fc,„ = Ll km (I 2m © BFT k )X^ m L\™{l m © Z^ 1 © I 2 (k-i)m)Dm m 

•(DCT-l ro+1 © DST-l m _x © 4_! © (DCT-2 m © DST-2 m ))s ? 2 r f m , 

where B^ m and X™ m are defined in (fOj) and (fA~5|) . and 
2fe ™-' --diagfl/ S in^±- 1)7rS 

V TO 



D 



'm+l 



diag ( 1 / cos ^ 



2m/ o<j 



0<j<m-l 

dia. 



<m 



iag (l/ 



. (j + l)7T 

sin 



2m J o<j<m 



(A.6) 



Appendix B. Proof of Theorem 16.21 

Consider .4 = = C[a;]/2T 2fem (a;) with basis {V (x), Vi(x), . . . , V 2 km-i(x)}. 
The corresponding polynomial transform is 

diag (l/ cos (fc t/ /2 ^ ) ' D CT-4 2fcm = DCT-4' 2km . (B.l) 

By Theorem 13. 2[ the polynomial r(ir) = T 2 k{x) generates the subalgebra 
B = (r(x)) = C[y]/2T m (y). 
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By Theorem l4.51 the £?-module VAx)B = C[y] / 2T m (y) for any < j < k. If we choose 



the basis {V^y)} 



0<£<m' 



then the polynomial transform is 



(fc + l/2K 
T f (cos ) 



= DCT-3 r . 



0<k,e<m 



Similarly, the 6-module W 3 {x){V 2k -i{x) - V 2k {x))/2 ■ B = C[y]/2T m (y) for any 
< j < k. If we choose the basis If we choose the basis {Ue(y)} 0<e<m , then the 
polynomial transform is 

W fc + 1/2)7r )1 (* + V2V 



0<fc,£<m-l 



diag ( 1/ sin y " ' ) -DST-3 m = DST-3' 

V TO /0<fc<m 

(B.2) 

We can verify using Theorem 14.61 that t 2 j — Vj(x) and t^j+i — Wj(x)(V 2 k-i(x) — 
V 2 k(x))/2 for < j < k, is a transversal of £? in .4. Hence, by Theorem 15. 11 we obtain 
the decomposition 

DCT-4' 2fem = M (I k ® (DCT-3^ © DST-3JJ) B. 

Here, B is the base change matrix from \x e \„^„^ , to the concatenation of bases 
of tj(x)B, < j < 2k, and by construction 

'1 



B = h 



Ll {m ^ • I m -x (g) DFT 2 



2fcm\T 
2m J • 



1; 



M is constructed as follows. Let 



M 



Let £>,■ = diag (tj(eti)) for < j < 2k. Then 

V / 0<i<n 

M = (D M | D X M Q | D 2 M | • • • | D 2k ^M ) 
We can simplify matrix M. Let us introduce matrices 
/ Co / c 



4 C4) m = 



\so 



. 4 S4) w 



Cfc-l/ 



\-so 



"Sfe-l\ 



Cfc-l / 



(B.3) 



Here, q = cos — — 2r )(2t+i)n and sg — sin ^ — 2r)(2i+i)7r ^ ^pj^ggg ma t r ices are used for 
the so-called s/cew DCT and DST [32]. Further, let us define rW = (2i + l)/(4m) and 

' r (i) +2j 

- , if j is odd 



if j is even 

2- r ( ' ) +2j 
fc 



for < j < LfJ- In case k is odd, we also define r k _ t 
define diagonal matrices 



r C*)_i 



1. Finally, let us 



Df 4) (r W ) = diag (l/cos(rf 



0<j<k 



D (S 4 )( r «) = diag ( sin {2krf*)l cos (rf V/2) 



o<j<fc 
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Then M = K%ML% m , where M 



( eS^f 4) (»" (i) )DCT-4 fc (rW)xf 4) (rW) e™" 1 ^ V l) )DST-4 fe (r«)4 S V ) ' 
^^^(^JDOT^frWjX^frW) 0^ 1 Di W) (rW)DST^(rW)4 S4) (rW) y 

We can further simplify ([B7T]| by substituting DCT-4 2fem and DST-3^ n with DCT-4 2fc „ 
DST-3 m using ()B.1|) and (|B.2I) . Then we use the equalities 



DST-3 m = diag [{-iy) -DCT-3 m -J m , 

\ / 0<j'<m 

DST-4* = diag f -DCT-Ak-Jk, 

V / 0<i<k 



to obtain the decomposition 



where 



DCT-4 2fem = Kl km {Kl m <g> DCT-4 fe )y^ m (DCT-3 m ® Lf){Kf k m ) T 
■h ® f if - 1} ■ Z m _! ® DFT 2 ] (^ m ) T , 



m ~ S U C4) (1 - r«)) (-l)^ 1 • J, • X<°*>(rW>) J (R4) 



and x£ C4) (r) is defined in (lB~3l) 
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